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A. 


One-Dimensional  and  MHD  Code  Techniques  for  Imploding  Plasmas  with 


Radiative  Energy  Flow 
1.  Motivation 

Existing  models  of  imploding  wire  or  puff  gas  plasmas^  leave  a 
number  of  physical  processes  unexplored  or  only  partially  examined  in  their 
assessment  of  diode  Imploded  plasmas  as  radiation  sources.  The  one-dlmension- 
al  Implosion  code  to  be  described  here  Is  an  attempt  to  Increase  the  physical 
relevance  of  such  models  and  to  enhance  the  versatility  and  robustness  of  the 
numerical  schemes  used  to  implement  them. 

Our  MHD  computational  efforts  are  aimed  at  improving  the  treatment  of  a 
number  of  processes  over  the  MHD  codes  presently  available.  These  include 

•  Photon-plasma  energy  exchange. 

•  Plasma  thermoelectric  radial  stresses  -  how  is  the  interior  of  the 
annular  plasma  pushed  by  the  J  x  B  forces  on  the  outer  surface? 

•  Current  and  field  duffusion  into  the  plasma  "skin"  -  is  the  current 
penetration  anomalously  rapid  as  in  classical  pinch  experiments, 
and  how  does  the  current  layer  evolve,  quantitatively.  The  details 
affect  the  extent  to  which  non-classical  corona  currents  reduce  the 
force  on  the  bulk  of  the  plasma. 

•  The  fractions  of  electromagnetic  power  absorbed  and  reflected  by 

the  moving  plasma  load  in  the  diode  cavity.  This  is  a  central  question 
since  it  measures  the  coupling  efficiency  of  the  hydrodynamic  driver. 
Can  the  experimental  diode  waveforms  be  correlated  to  the  plasma 
motion? 

•  The  dynamic  energy  partition  between  electron  thermal  energy  and 
ionization  potential . 
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These  considerations  can  be  organized  into  several  concrete  goals  for  the 
calculation. 


1)  Find  how  important  are  the  details  of  photon/plasma  energy 
exchange  in  structuring  the  implosion  trajectory  and  (spatially 
resolved)  energetics.  The  transitions  in  qualitative  characteris¬ 
tics  (ranging  from  refrlgerative  collapse  to  strong  reflection) 
noted  in  the  core-corona  model®  were  quite  sensitive  to  the  opacity 
calculation.  In  a  full  1-D  implosion  calculation  the  resolution 

of  the  enhanced  "thermal  conduction"  due  to  direct  photon/plasma 
energy  exchange  will  provide  a  more  accurate  picture  of  the  plasma 
emission  profile. 

2)  Determine  the  role  of  separate  ion  and  electron  temperatures 
in  establishing  this  emission  profile  and  in  governing  the  plasma 
dynamics.  As  compared  with  a  single-temperature  picture,  one 
expects  the  proper  two-temperature  treatment  to  alter  the  pressure 
gradients  in  low  density  regions  and  change  the  opacity  in  higher 
density  regions,  perhaps  altering  the  spectral  mix  of  line  and 
continuum  radiation. 

3)  Examine  the  effect  of  the  chemical  potential  (properly  included 
in  the  plasma  energy  balance  as  a  function  of  electron  energy  and 
ion  density)  on  the  implosion  energetics  and  the  emission  profile.4 
An  accurate  partition  of  incoming  energy  will  probably  provide  lower 
(more  physical)  temperatures  at  peak  implosion  than  earlier  models. 

4)  Develop  an  accurate  picture  of  the  current  penetration  mechan¬ 
isms.  Diffusive  transport  of  E2(r ,  t)  (or  )  is  dominant  to  a  first 
approximation,  but  overall  implosion  performance  is  sensitive  to  the 
effective  skin  depth.®  The  role  of  the  free  wave  field  components 


in  setting  this  skin  depth  has  not  been  examined,  and  the  rate  of 
Inward  propagation  for  the  drift-speed  limitations  In  J2  maty  be 
affected  by  these  components. 

5)  Resolve  the  sources  of  reflected  power  within  the  plasma  load 
and  relate  them  to  the  available  measured  waveforms  from  experiment. 

An  electromagnetic  calculation  of  all  fields  within  the  plane  parallel 
diode  opens  one  more  channel  of  diagnostic  Information,  and  may 
provide  a  means  of  inferring  plasma  motion  from  measured  waveforms. 

It  will  also  provide  a  more  accurate  statement  of  energy  conservation 
because  the  energy  contained  in  the  net  free  wave  component  will  be 
included. 

6)  Examine  the  effect  of  radial  ambipolar  electric  fields,  and  the 
resulting  stress,  on  the  implosion  dynamics,  particularly  in  the  low 
density  regions  where  drift  speed  limits  are  expected.  The  radial 
electric  fields  arise  from  thermoelectric  effects  and  axial  drifts 
of  charge  carriers.  Large  temperature  gradients,  allowed  by  weak 
cross  field  thermal  conduction  and  enhanced  ohmic  heating  in  the 
corona  plasma,  may  provide  a  significant  radial  electric  field  and 
resulting  inward  stress  on  the  plasma. 

Once  these  physical  questions  are  addressed  coherently,  a  smooth 
upgrading  in  the  sophistication  of  the  radiation  emission  and  transport 
model  becomes  a  fruitful  exercise.  With  the  radiative,  hydrodynamic,  and 
electrodynamic  models  structured  from  first  principles  there  are  no  "free 
knobs"  and  the  subsequent  benchmarks  with  an  appropriately  configured 
experiment,  if  successful,  provide  a  firm  basis  for  similar  theory  in  more 
difficult  and  perhaps  more  practical  diode  configurations. 
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Apart  from  these  physical  questions  are  various  topics  arising  In 
the  numerical  implementation  of  whatever  implosion  model  Is  selected.  First 
is  the  choice  of  grid;  a  typical  simulation  will  show  the  model  plasma 
compressed  several  orders  of  magnitude  during  the  course  of  a  calculation. 

A  strict  Eulerian  fluid  calculation  can  thus  lose  its  spatial  resolution 
at  the  final  collapse,  when  this  resolution  is  most  needed.  Moreover  large 
density  gradients  are  conmionpl ace ,  and  it  is  probable  that  the  dynamics  of  the 
low  density  region  are  dominated  by  the  balance  of  density  gradients  and 
magnetostri cti ve  stresses.  A  Lagrangian  description  is  thus  preferable  in 
an  implosion/explosion  calculation,  so  long  as  the  external  stresses  on  the 
system  can  be  defined  readily  over  the  spatial  domain  of  interest.  On  the 
other  hand  the  equations  governing  the  field  (Ez)  generating  the  implosion  are 
solved  most  naturally  in  the  fixed  laboratory  frame.  The  electromagnetic 
problem  therefore  requires  an  interpolation  capability  which  can  map  external 
stresses  onto  a  convecting  fluid  and  extract  field  sources  from  this  moving 
fluid.  The  natural  compromise  adopted  here  is  a  completely  Lagrangian  fluid 
model  convecting  through  a  laboratory-frame  (Eulerian)  electromagnetic  grid. 
The  electromagnetic  grid  can  be  made  sufficiently  non-uniform  to  concentrate 
the  Ez  information  density  near  the  axis  while  the  smoothly  distorted 
Lagrangian  fluid  grid  is  a  natural  choice  for  an  accurate  resolution  of  density 
gradients. 

The  second  numerical  topic  is  the  choice  of  Ez  or  Bg  as  the  fundamental 
el ectrodynamic  variable.  In  our  view  the  more  convenient  selection  is  Ez 
because  it  is  the  natural  field  from  which  to  establish  drift-speed  limitations 
on  the  current  density.  Using  Bg  and  J2  =  ^  (V  *  Bq)2  in  a  drift  speed 
limitation  algorithm  has  two  disadvantages:  the  current  is  derived  from  a 
potentially  noisy  differentiation  process,  and  the  limit  J2  =  eng  c$  does  not 
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imply  a  completely  local  value  for  Ez  because  the  conductivity  depends  on 
Bfl,  and  hence  on  J,  at  other  grid  points.  It  is  clearly  possible  to  remove 
the  second  problem  by  iteration,  but  the  physical  process  providing  the  drift 
speed  limit  Is  a  local  (fine  scale)  one.  Proceeding  first  from  the  local 
values  of  Ez,  then  iteratively  establishing  small  corrections  to  the  (possibly 
limited)  Jz  from  non-local  considerations,  is  a  more  natural  way  to  Impose 
the  limit  and  is  probably  less  noisy. 

A  third  point  is  the  selection  of  forward  time  integration  methods. 

A  coimon  preference  for  conservatively  advanced  explicit  integrators  must 
be  examined  carefully  in  the  context  of  this  radiatively  coupled  problem. 

As  the  level  of  sophistication  in  the  radiative  transfer  physics  is 
Increased  a  considerable  expense  is  incurred  for  each  fluid  time  deriva¬ 
tive.  This  expense  will  scale  linearly  with  the  number  of  frequency  domains 
and  quadratically  with  the  spatial  resolution.  In  order  to  optimize  the 
complete  algorithm  it  is  therefore  necessary  to  make  the  best  possible  use 
of  every  time  derivative  available  from  the  space-time  mesh.  Moreover,  as 
distinct  from  calculations  dominated  by  processes  internal  to  the  fluid,  the 
problem  of  simulating  an  electromagnetically  driven  implosion  requires  a 
self-consistent  time  evolution  of  both  field  and  fluid.  This  class  of 
problems  is  more  likely  to  be  handled  efficiently  and  accurately  by  an 
implicit  method,  or  at  least  by  a  predictor/corrector  scheme.  This  choice 
has  only  the  disadvantage  of  being  more  cumbersome  to  implement,  and,  as 
described  below,  the  methods  chosen  make  a  smooth  transition  from  "nearly 
explicit"  to  "fully  implicit."  As  dictated  by  the  performance  characteris¬ 
tics  of  the  complete  algorithm,  the  user  will  be  able  to  optimize  the 
accuracy  and  cost  tradeoff  by  several  means. 
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A  final  consideration  Is  the  accuracy  of  spatial  differentiation 
within  either  the  hydrodynamic  or  electrodynamic  algorithms.  If  these 
algorithms  can  make  good  use  of  non-uniform  meshes  in  space  and  time,  then 
one  can  sustain  useful  accuracy  with  fewer  points.  The  common  differentiation 
methods  based  on  finite  differences  are  quite  accurate  when  a  dense  mesh  is 
adnissible,  but  if  one  wishes  to  keep  down  the  cost  of  radiative  transfer 
calculations  a  higher  order,  smoother  scheme  is  preferable.  Some  useful 
alternatives  are  discussed  in  Section  A.3.b.  We  have  developed  such  schemes, 

involving  smooth  interpolators ,  as  one  segment  of  our  computational  effort, 
and  this  has  already  shown  its  utility  in  interfacing  various  parts  of  the 
NRL  radiation  physics  with  our  plasmadynamic  codes. 


2.  ID  MHD  Code:  Computational  Approach 

The  specific  physical  content  of  the  model  has  been  described 
6  7 

elsewhere.  *  The  task  of  implementation  can  be  readily  broken  down  into 
the  three  problem  areas  treated  below  -  electromagnetics,  hydrodynamics, 
and  radiative  transfer.  The  requirements  for  each  problem  area  and  a 
summary  of  the  algorithms  employed  are  given. 

The  electrodynamic  problem  statement  is  simple:  given  voltage  wave¬ 
forms  at  the  periphery  of  the  diode,  establish  the  E  and  B  fields  everywhere 
in  the  interior  for  the  duration  of  the  implosion  process.  By 
driving  the  diode  E/M  model  at  one  radius  with  a  voltage  waveform,  the  current 
drawn  by  the  load  and  the  voltage  waveforms  at  any  interior  point  are  avail¬ 
able  as  predicted  observables/diagnostics  of  the  code  and  the  experiment. 
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Conventional  approaches  replace  the  vacuum  portion  of  the  diode  with  a 
circuit  relation  derived  from  Faraday's  law,  but  a  complete  electromagnetic 
calculation  provides  more  information.  If  this  latter  goal  is  chosen  then  the 
issue  of  efficient  and  effective  calculation  becomes  an  important  question. 
With  the  implosion  model  discussed  here  the  choice  has  been  to  admit  either 
option,  the  full  electromagnetic  computation  or  the  diffusive  limit  and 
circuit  equation. 

For  the  first  option,  it  is  instructive  and  in  many  ways  practical 
to  utilize  the  generalized  Hertz  vector  potential.  First,  this  new  potential 
reduces  the  electrodynamics  to  a  single  component  of  a  single  vector  wave 
equation.  This  reduction  simplifies  the  calculation  by  placing  no  self- 
consistency  constraints  on  the  time  evolution  of  the  fields  E  ,  BQ.  The 
required  self-consistency  is  automatic  when  E„,  Bfl  are  derived  from  a 
single  potential  Z2-  Second,  the  coupling  term  between  field  and  plasma 
is  properly  defined  to  all  orders  in  the  radial  fluid  convection  velocity 
and  smoothly  admits  marginally  stable  drift-speed  limited 
conduction  and  (quasi-static)  thermo-electric  contributions.  Properly,  as 
fast  time  scale  phenomena  which  track  the  instantaneous  state  of  the  system, 
these  effects  should  be  calculated  from  a  single  time  ?lice  of  information. 
Without  the  Hertz  vector  formulation,  an  electromagnetic  calculation  evolving 
Ez'  Be  seP»r«tely  can  only  accomplish  this  to  the  extent  that  it  satisfies 
the  temporal  self-consistency  constraint  mentioned  above.  Third,  the  coupling 
term  is  a  spatial  convolution  with  a  fixed  Green's  function.  The  spatial 
averaging  will  generally  serve  to  lower  the  noise  level  of  the  calculation 
by  damping  any  chatter  that  may  appear  on  the  discrete  source  array  extracted 
from  the  plasma  representation.  Fourth,  while  the  use  of  Z  does  reouire 
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extraction  of  high  order  spatial  derivatives  from  discrete  data  to  form  the 
E  and  B  fields  and  the  sources  for  the  Z  equation,  the  accuracy  of  such  spatial 
differencing  can  be  improved  by  the  use  of  so-called  smooth  interpolators.8 
The  estimation  of  local  Ez  fields  from  the  parent  Hertz  vector  by  means  of 

these  splines  is  routinely  accurate  to  [10’6],  even  between  the  input  data 
poiqts. 

The  Hertz-vector-based  electromagnetic  calculation  proceeds  from  a 
spatially  discrete  set  of  information  {Z,  3TZ,  J,  3TJ,  3^J},  given  at  a  single 
time.  The  first  two  elements  of  this  set  are  defined  at  fixed  radial  points; 
the  last  three,  established  (sequentially)  on  a  convecting  (Lagrangian)  grid  of 
plasma  cell  centers.  The  local  values  of  E,  B  at  any  such  time  level  are  deter¬ 
mined  by  (Z,  3tZ,  J)  and  partially  determine  the  plasma  stress  and  heating  rates. 

2 

From  the  values  of  3t«J,  a  spatial  convolution  produces  the  Hertz  vector 

source  terms.  These  source  terms  are  interpolated  for  use  in  an  explicit 
wave  quadrature  formula  which  advances  the  Hertz  vector  and  its  partial  time 
derivative  in  order  to  establish  Z,  S^Z  at  any  subcycle  time  step  (used  in 
the  thermal  conduction  scheme)  or  at  the  next  major  time  level. 

The  second  electrodynamic  option  is  that  of  E-field  diffusion  (which 
can  be  derived  from  the  Hertz  potential  wave  equation  in  the  limit  that  incom¬ 
ing  and  outgoing  waves  are  in  detailed  balance).  This  limit  is  approximate  but 
reasonably  accurate  in  the  plasma  where  the  conduction  current  exceeds  the 
displacement  current.  It  neglects  the  net  displacement  current,  making  any 

detailed  outgoing  wave  calculation  impossible  because  no  wave  source  is 

extracted  from  the  plasma.  Reflected  power  is  manifested  solely  as  a 

diminished  voltage  on  the  plasma  load  and  the  energy  delivered  is  simply 

/^dt  (I  ,  *V  ).  Here  only  the  plasma  current  waveform  is  available  as 
plasma  gen 
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diagnostic  Information,  Insofar  as  the  voltage  at  the  convectlng  plasma 
vacuum  Interface  Is  not  usually  available  f .  m  experiment.  The  major  simpli¬ 
fication  Is  that  Z  need  not  be  explicitly  employed;  It  Is  sufficient  to 
diffuse  the  E„  field  and  calculate  the  B.  field  in  a  (.self-consistent)  quasi- 
static  manner.  Diffusion  of  E,  rather  than  B.  Is  preferable  because  the 

Z  w 

criteria  for  drift  speed  limitation  stem  most  directly  from  the  application  of 
an  "E-field-to-local-O"  mapping. 

Diffusion-based  electrodynamics  proceeds  concurrently  with  subcycled 

thermal  conduction.  A  circuit  relation  sets  the  time  varying  value  of  Ez  at 

the  outer  boundary  of  the  plasma,  and  the  time  slice  information  {E^,  E',  B, 

8,  a,  a}  defined  on  the  (convecting)  plasma  cell  centers  determines  the 

material  derivative  of  E  .  Here  E’  3  E  +  8„  x  Ba  is  the  field  in  the 

z  z  r  0 

convecting  frame;  Eth  is  the  axial  thermoelectric  field,  arising  from  radial 
pressure  gradients.  The  new  field  E2  at  the  next  time  level  is  determined 
implicitly  by  the  time  sequence  of  material  derivatives.  The  plasma  mesh 
evolves  in  response  to  the  stresses  partially  determined  by  the  local  values 
of  E  and  B.  In  this  mode  of  operation  some  iteration  is  required  in  obtaining 
both  DEZ/Dt  and  BQ,  in  contrast  to  the  electromagnetic  calculation.  (B0 
requires  iterative  refinement  due  to  the  magnetic  field  dependence  of  the 
plasma  conductivity.) 

In  either  the  electromagnetic  or  the  diffusive  mode,  the  imposition 
of  drift-speed  limited  conduction  is  straightforward  because  the  E2(r,t)  pro¬ 
file  is  established  a  priori,  before  J.  is  estimated.  The  rational  basis  for 
preference  between  the  two  modes  can  be  derived  only  when  the  goals  of  item  (4) 
above  are  reached  and  when  the  comparative  computation  costs  of  well -tuned 
subroutines  are  established. 
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The  hydrodynair.ic  problem  Is  readily  formulated  (once  E,  B  are  con¬ 
sidered  given)  In  terms  of  well-known  techniques.  The  primary  challenge  is 
to  produce  a  quiet,  relatively  robust  algorithm  capable  of  concentrating  the 
computational  effort  on  the  faster  time  scales  of  the  problem,  minimizing  the 
number  of  expensive  derivative  evaluations  required,  and  retaining  spatial 
resolution  of  very  highly  compressed  plasma  loads. 

In  view  of  these  requirements  a  Lagrangian  formulation  of  the  two- 
fluid  transport  equations  for  number  density,  flow  velocity,  and  temperature 
has  been  implemented.  There  are  several  simplifications  to  the  complete  set 

of  transport  equations  which  are  appropriate  to  the  imploding  diode  plasma, 

9 

but  the  hydrodynamic  model  is  essentially  that  of  Braginskii.  In  particu¬ 
lar  the  radial  flow  field  is  assumed  to  be  identical  for  ions  and  electrons, 
the  electron  density  is  assumed  to  be  in  quasineutral  equilibrium  with  the  ion 
density  (allowing  only  small  ambipolar  charge  separations),  and  the  system  of 
moment  equations  is  closed  with  an  ideal  gas  equation  of  state  relating  the 
pressure  and  temperature  for  each  fluid.  This  reduces  the  fluid  variables  to 
(rij,  V^,  Tj,  Tg»  £j)where^is  the  degree  of  ionization  and  the  average 
ionization  potential  for  the  current  Rvalue.  However  Tg  is  not  a  very  good 

thermodynamic  variable  because  it  is  tightly  coupled  to  and  to  £j.  In 

2  £I 

order  to  follow  the  ionization  dynamics  efficiently,  we  let  9  *  T  +  y-s — 

e  e 

and  evolve  9g  instead.  The  approximation  of  the  collisional  radiative 
equilibrium  (CRE)  model^  is  to  assume  a  very  rapid  establishment  of  ^  and 
from  (given)  nj  and  T  .  Hence  the  usual  heating  sources  and  sinks  for  T  , 
when  ^  is  fixed,  become  (in  effect)  heating  sources  and  sinks  for  9g  in  the 
context  of  CRE  ionization  dynamics.  From  a  sufficiently  accurate  representa¬ 
tion  °f  ^(0e)  and  £j(9g)  preprocessed  by  the  complete  CRE  model  and  smoothly  inter- 
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polated,  one  may  simply  advance  6  Cthe  "grand  canonical  temperature")  and 
compute  Tg,  ^  ,  Cj  and  all  other  transport  coefficients  from  the  0-value 
obtained  at  any  timestep.  This  latter  method  is  presently  implemented.  As  is 
common  practice  in  Lagrangian  calculations,  one  may  also  remove  the  ion 
density  from  the  set  of  evolved  quantities  by  conserving  the  total  number  of 
ions  in  any  cell,  unless  a  rezonfhg  is  done,  and  calculating  the  density  n^ 
from  the  time  varying  cell  volume.  Hence  the  irreducible  set  of  fluid  vari¬ 
ables  to  be  evolved  is  {V  ,  Tt,  6  },  with  the  r  5  dt.Vft.)  defining  the 
variable  mesh  positions  to  be  used  in  computing  the  density  nj.  For  the 
choice  of  diffusive  electrodynamics,  this  set  is  expanded  to  {E  ,  V  ,  Tt  0  } 
and  the  separate  evolution  of  {Z,  Z}  through  the  Hertz  vector  wave  equation 
is  discarded. 

The  complete  set  of  magnetized,  ^-dependent  transport  coefficients 
has  been  coded  in  a  single  subroutine,  and  those  required  for  the  simplest 
relevant  problem  are  the  axial  electrical  conductivity,  the  radial  ion  and 
electron  thermal  conductivities  and  the  radial  thermoelectric  coefficient. 

The  expansion  of  the  model  to  include  radial  ambipolar  electric  fields, 
viscous  stresses  and  viscous  heating  only  requires  activation  of  the  axial 
thermoelectric  and  friction  coefficients,  and  five  viscosity  coefficients.  More¬ 
over,  because  the  relaxation  time  is  the  cornerstone  of  the  transport  calculation 
an  extension  into  more  strongly  coupled  plasma  domains  is  quite  straightforward 
on  a  cell -by-cell  basis.  The  use  of  a  modified  relaxation  time  based  on  the 
local  Coulomb  logarithm  is  a  useful  first  approximation  in  strongly  coupled 
plasma  transport^  and  can  be  relevant  here  in  some  simulations  of  wire  or  foil 
annular  implosions. 
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The  hydrodynamic  calculation  across  a  major  time  step  proceeds  from 
the  {Nj,  Vr,  Tj,  0e>  at  the  beginning  of  the  time  step.  Here  Nj  (.the  time  Invar¬ 
iant  total  number  of  cell  Ions),  Tj  and  0g  are  defined  at  the  cell  centers,  while 

V  and  r(t)  are  defined  at  the  cell  Interfaces.  From  this  basic  data  the  model 
r 

calculates  Te»  £j»  *11  transport  coefficients,  all  heating  rates  and  all 
fluid  stress  terms  (after  E  and  B  have  been  established).  The  major  time 
step  is  selected  from  a  competition  among  the  magnetosonlc  Courant  condition, 
fractional  cell  area  changes,  local  truncation  error  in  estimating  future 
accelerations  (DV^Dt),  and  the  most  rapidly  varying  electromagnetic  source 
terms.  Once  the  major  time  step  is  set  a  time  window  Is  defined  over  which 
thermal  conduction,  compressional  heating,  radiative  heati.jg  or  cooling,  and 
ion-electron  exchange  heating  are  subcycled  at  the  cell  centers  (under  an 
assumed  cell  boundary  evolution).  This  subcycle  Is  implemented  using  an 
(implicit)  variable  order  Adams  predictor-corrector  or  Gear  method.  At  the 
completion  of  the  subcycle  the  explicit  flow  field  advance  (and  implied  cell 
boundary  motion)  is  either  retained  or  iteratively  refined  (becoming  the 
starting  point  for  an  implicit  method).  If  desired  the  iteration  proceeds 
until  self-consistency  is  achieved.  The  self-consistency  criterion  is  solution 
of  the  first  order  (non-linear)  difference  equation  for  the  fluid  accelera¬ 
tion,  to  a  specified  tolerance  at  all  mesh  points.  For  the  diffusion  based 
electrodynamics  the  subcycle  includes  the  E  field  evolution  as  well. 

A  central  element  in  the  hydrodynamic  model  is  the  question  of 
radiation  transport.  The  algorithm  presently  implemented  is  a  compromise 
between  the  full  CRE  calculation  and  the  much  simplified  local  approxima- 
tion  of  the  SPLAT  code.  First,  the  emission  function  preprocessed  by  the 
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CRE  model  is  represented  over  the  relevant  density  and  temperature  domain. 

One  such  emission  curve  Is  required  for  each  radiation  cateogry  one  wishes 

to  calculate,  A  simple  set  of  radiation  categories  Is:  {lines  hv  <  1  keV, 

lines  hv  >  1  keV,  free-bound  continuum,  &  Bremsstrahlung).  Next,  for  each 

13 

radiation  category  a  damping  coefficient  Is  computed  for  each  cell.  These 
damping  coefficients  are  path  Integrated  along  a  selected  ray  to  produce  a 
probability  of  escape  for  that  photon  category.  From  the  escape  probability 
field  a  matrix  of  coupling  coefficients  is  derived.  Once  the  coupling  matrix 
is  given,  selected  inner  products  with  the  vector  of  cell  emissions  for  a 
particular  category  produce  the  cell-to-cell  photon  exchange  and  the  net 
radiative  loss  to  the  plasma  from  any  zone  in  each  radiation  category. 

This  formulation  has  the  advantages  of  being  easily  expanded  into  a 
multiple  group  algorithm  and  easily  modified  to  include  other  elements.  At 
present  At  and  Ar  are  available,  although  the  representation  of  9e  is  not 
yet  available  for  Ar.  The  disadvantages  lie  in  the  added  calculation  of  the 
coupling  matrix  and  in  the  need  to  establish  a  fresh  matrix  for  each  radiation 
category.  It  is  hoped  that  the  use  of  implicit  methods  for  the  overall  hydro- 
dynamic  advance  will  offset  the  cost  of  these  radiative  loss  calculations  by 
requiring  fewer  time  derivative  evaluations  to  achieve  useful  accuracy. 
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3.  Summary  of  Code  Modules  Developed 


a.  The  Imploding  Plasma  Radiator  Code  (Program  ZDIPR) 

A  computer  program  Is  being  constructed  to  solve  the  one-dimensional 
(in  radius)  time-dependent  interaction  of  a  plasma  annulus  (with  temperature- 
dependent  degree  of  ionization)  and  an  electromagnetic  driving  field  In  a  cy¬ 
lindrical  cavity.  It  was  shown  previously  that  the  number  of  equations  describ¬ 
ing  the  fields  can  be  reduced  and  the  computational  stability  probably  improved 
by  the  use  of  a  generalized  Hertz  vector,  Z,  from  which  E  and  B  field  vectors 
can  be  derived  by  differentiations.  In  the  limit  of  high  plasma  density  the 
wavelike  term  in  the  field  equations  can  be  neglected;  this  corresponds  to 
a  field-diffusion  approximation.  In  simple  cases  —  constant  scalar  con¬ 
ductivity  for  example— this  gives  the  usual  magnetic  diffusion  equation. 

The  program  must  couple  field  or  Hertz  potential  quantities, 
evaluated  most  conveniently  on  a  stationary  (Eulerian)  grid  of  points,  with 
fluid  equations  describing  the  mass-motion  (v)  electron  and  ion  temperatures 
(Tfi,  Ij)»  ionization  (j),  and  effective  conductivity  (cr)  of  the  plasma,  all 
functions  of  r  and  t  evaluated  most  conveniently  on  a  Lagrangian,  or  co¬ 
moving  grid.  The  differential  equations  are  integrated  forward  in  time  by  a 
variable-timestep  method  (GEAR),  but  because  the  plasma  thermodynamic  state 
evolves  much  more  rapidly  than  the  fields  do,  it  is  updated  separately,  on  a 
fast  timescale,  with  subcycle  timesteps  shorter  than  the  major  timestep  T 
required  for  integration  of  the  field  and  momentum  equations.  Doing  these 
updates  on  different  timescales  is  expected  to  save  a  significant  amount  of 
running  time  and  expense  in  the  computation,  and  may  make  the  difference 
between  a  practical  cost  and  an  impractically  large  one. 
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In  the  code  the  subroutines  that  update  the  plasma  state  on  the 
Lagranglan  "r"  mesh  are  MESHSTRESS  (forces),  TOOT  (heating),  CREOS  (Ioniza¬ 
tion  state)  and  FLUIDOTS  (material  derivatives),  all  of  these  being  used  by 
the  subroutine  GEARBOX  which  controls  the  updating  as  part  of  the  overall 
fluid-advance  routine  HYOROPUSH. 

The  subroutine  that  updates  the  field  or  Hertz  potential  variables  on 
the  Eulerlan  "RM  mesh  Is  called  EPUSH  or  ZPUSH.  When  the  Hertz  potential  repre 
sentatlon  Is  used,  the  subroutine  HERTZDOTS  updates  the  source  terms  In  the 

potential  equation  at  each  major  timestep.  When  the  electric  field  diffusion 
representation  is  used,  the  corresponding  subroutine  EDOT  updates  the  field 
source  $TJ  and  alternate  plasma  state  movers  (.GEARBOX  +  TETGEAR;  FLUIDOTS  ■* 
TETDOTS)  are  used  because  different  electromagnetic  information  Is  required. 

This  code,  called  ZDIPR,  is  now  mostly  written  and  some  of  its 
portions  have  been  separately  checked.  Further  testing  and  integration 
of  the  various  subroutines  will  be  carried  on  in  the  coming  weeks. 


b.  The  Smooth  Interpolator  Packages 

Computer  programs  have  been  written  (and  are  currently  in  use) 
which  fit  smooth  spline  functions  to  discrete  data.  The  smooth  functions 
are  sums  of  Gaussians,  error  functions,  or  integrals  of  the  error  function, 
at  the  user's  option.  The  user  specifies  a  typical  width  parameter,  and  the 
smoothness  of  the  representing  functions  depends  on  the  width  parameter  chosen 
(iteration  of  the  choice  of  width  may  be  necessary).  When  optimized,  the 
fitting  functions  also  have  smooth  first  and  second  derivatives  (to  within  6 
significant  figures  in  a  typical  example).  These  programs  -  REPIT  which 
represents  the  data  by  spline  functions  (using  subroutine  SMILAM),  and 
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CHECKIT,  which  evaluates  the  interpolant  over  any  specified  domain,  are 
described  more  fully,  along  with  the  underlying  theory,  In  Chapter  II  of 
the  accompanying  Technical  Report  (.JAYCOR  Tech.  Report  J207-82-002) .  These 
user-oriented  programs  are  of  great  Importance  In  the  plasma-radiation  and 
plasma  dynamics  community,  and  have  been  made  available  to  users  working  in 
these  areas  at  NRL.  They  are  a  necessity  In  codes  which  are  required  to 
extract  smooth  fields  from  potentials  given  at  discrete  data  points. 
Including  the  ZDIPR  code  described  in  Section  3. a.  above. 
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4.  Status  of  ID  MHD  Code 


The  HYDROPUSH  ensemble  of  subroutines  has  been  partially  written; 
given  a  predetermined  E-fleld,  the  time  step  choice  STEPPER  and  the  iteration 
of  the  mesh  are  coded.  FLUIDOTS,  which  determines  the  material  derivatives 
of  velocity  and  temperature,  has  been  coded;  its  companion  for  the  diffusive 
code  has  been  designed  and  is  being  coded.  GEARBOX,  which  carries  out  the 
integration,  has  been  coded.  The  Jacobian  routines,  which  assist  DGEAR  in 
integrating  the  equations  efficiently,  are  being  designed. 

The  field  solver  ZPUSH  has  been  coded  and  the  capability  of  unraveling 
E  and  B  from  the  generalized  Hertz  potential  Z  has  been  established.  (EFIELD 
and  BFIELD  are  the  subroutines  which  do  this.)  The  evolution  of  Z  in  time 
is  accomplished  by  WAVEQUAD  and  six  functions  which  determine  the  elements 
of  the  wave  equation  quadrature's  intergrand.  The  code  WAVEQUAO  has  been 
partially  written  and  is  undergoing  testing  at  this  time.  The  coding  of 
EDOTS,  which  couples  the  external  circuit  to  the  plasma  "boundary"  current, 
is  beginning.  These  subroutines  all  rely  on  previously-developed  and  tested 
routines,  e.g.  CREOS,  BRAGINSKY,  TOOT,  which  evaluate  the  basic  plasma  and 
radiation  quantities. 

A  broad  spectrum  of  smooth  interpolator  routines,  just  described  in 
Section  A.3.b.  and  the  accompanying  technical  report,  has  been  built, 
reduced  to  practice,  and  put  to  use  in  the  EFIELD  and  BFIELD  routines. 
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B.  Simple  Models  for  Radiation.  Dynamics,  and  Generator  Optimization 

The  search  for  simple  but  approximate  models  capable  of  giving 
radiation  scaling  laws  at  low  cost  continues  alongside  the  larger-scale 
computer  modeling.  To  this  end  we  have  Initiated  three  interrelated  calcula¬ 
tions: 

(1)  radiation  blanketing  and  energetics  for  cylinders  with  Gaussian  ion 
density  profiles. 

(2)  generalizing  the  JAYCOR  SQUEEZE  code,  which  describes  approxi¬ 
mately  the  dynamics  of  the  plasma  and  radiation  after  assembly  onto 
the  axis,  and  interlinking  it  with  the  SAI  WIRES  code,  which 
describes  primarily  the  run-in  phase  of  the  implosion  prior  to 
assembly.  The  results  of  calculation  Cl)  are  to  be  used  in  the 
equations. 

(3)  a  calculus-of-variations  approach  to  deriving  explicit  equations 
or  computational  algorithms  for  the  machine  current  or  voltage 
pulse  profiles  and  the  plasma  initial  conditions  that  optimize  the 
desired  radiative  loss.  The  results  of  calculations  Cl)  and  C2)  are 
to  be  used  in  this  one. 

Each  of  these  areas  will  be  described  briefly  now  to  show  their  justification 
and  interrelation  as  well  as  the  mode  of  approach  to  them. 

1.  Radiation  Functions  for  Gaussian  Cylinders 

An  assembled  radiating  plasma  cylinder,  partly  opaque  to  its  own 
radiation,  is  described  by  a  density  profile  n.(r)  with  typical  radius  R 
(say  the  radius  at  half-maximum  of  1/e-maximum  density),  a  central  density 
nQ  *  n^o),  and  an  electron  temperature  profile  T(r).  If  radiative  heat 
transport  were  very  fast,  one  could  imagine  T(r)  being  uniform:  T(r)  *  T 
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everywhere.  In  that  oversimplified  baseline  case  one  can  use  the  NRL  radiative 
transport  multicell  code  to  calculate  the  total  radiative  energy  loss  rates, 
for  photons  with  energy  above  and  below  energy  E  ,  i.e. 

P>  (n0,  tni ( r) ] ,  T),  ?<  CnQ,  T). 

Because  most  Interesting  density  profiles  qualitatively  appear 
Gaussian,  we  specialize  in  this  baseline  case  to  Gaussian  n^(r),  so  that  the 
radiative  loss  functions  become 

P*  ln0,  R,  T),  P<  (nQ,  R,  T). 

However,  in  reality,  at  the  edges  of  such  plasmas  the  radiative  heat  transport 
can  never  keep  the  plasma  at  the  central  temperature.  Radiative  heat  loss 
insures  that  there  will  be  a  fall  off  of  T (r),  i.e.  a  cooler  blanket  surrounding 
the  hot  region.  This  blanket  generates  radiation  in  the  P<  category,  and 
absorbs  radiation  in  the  P*  class,  so  it  is  important  in  both  the  diagnostics 
and  the  energetics.  Even  if  the  blanket  layer  is  thin  compared  with  the 
cylinder  radius,  its  optical  depth  may  be  significant  (as  in  the  case  of 
stars),  because  the  local  opacity  increases  as  the  temperature  falls  with  r. 
There  may  also  be  a  hotter  tenuous  corona  outside  the  blanket,  but  the  density 
is  too  low  to  be  of  importance  in  the  radiation  energetics. 

So  we  are  preparing  to  calculate  with  multicell  techniques  the 
corresponding  self-consistent  T(r)  and  P<  and  P*  for  a  Gaussian  density  profile 
(with  uniform  heat  input  maintaining  T(.r)  constant),  and  thus  to  characterize 
the  width  of  the  blanket  6R  (nQ,  R,  T  )  formed  by  the  temperature  profile 
T(r;  nQ,  R),  as  well  as  the  new  radiated  powers 
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P>  (nQ,  R,  T0)  and  P<  (nQ,  R,  TQ), 

where  Tq  Is  the  central  temperature. 

Simple  analytic  fits  to  these  functions  will  then  provide  a  useful 
tool  in  that  they  can  be  put  into  the  energy  equation,  which  together  with 
the  momentum  equation  and  a  circuit  equation  determines  the  radiation  scaling 
of  such  a  model  assembled  plasma.  The  momentum,  energy,  and  circuit  equations, 
in  such  a  simple  model,  are  three  ordinary  differential  equations,  and  their 
computational  solution  is  much  simpler  than  the  MHO  fluid  equations  used  for 
the  more  exact  modeling.  These  simple  equations  and  their  method  of  solution 
are  discussed  next. 

2.  Interlinked  WIRES  and  SQUEEZE  Codes 

In  a  combined  JAYCOR/SAI/NRL  effort  characteristic  of  this  program  we 

14 

have  begun  combining  the  SAI  WIRES  code  (Mondell i,  1981)  for  the  run-in 
phase  of  a  multi-wire  array  implosion,  with  a  generalization  of  the  JAYCOR 
SQUEEZE  code  (Guillory,  1980)  for  the  assembled  plasma  compression  and  radiation 
phase.  The  WIRES  code  model  is  followed  until  the  imploded,  thermally  expanded 
wires  just  touch.  The  plasma  is  then  smoothed  to  a  thick  annular  layer,  and 
continues  the  implosion  in  a  manner  akin  to  the  JAYCOR  SPLAT  code  until  a 
non-hollow  cylinder  is  formed.  The  SQUEEZE  code  then  describes  the  dynamics 
of  the  assembled  plasma,  at  least  until  several  growth  times  for  macroscopic 
fluid  instabilities,  whereupon  the  plasma  no  longer  follows  any  simple  one¬ 
dimensional  theory. 

14 

The  WIRES  code,  described  more  fully  elsewhere,  consists 
of  a  simple  Newton's-law  dynamical  equation  for  the  radial  position 
of  the  N  wires,  a  circuit  equation  relating  the  voltage,  current,  and  time- 
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varying  inductance  of  the  array,  and  an  equation  for  heating  of  the  wires  in 
the  presence  of  blackbody  radiative  losses,  with  the  size  of  individual  wires 
swelling  as  they  heat.  When  these  cylinders  touch,  they  are  replaced  by  an 
annulus  of  plasma  of  the  same  total  area,  and  this  is  further  compressed  onto 
the  axis  by  the  total  current  in  such  a  way  that  the  inner  radius  b  moves  at 
the  finite  speed  |b  |  *  |a|  +  c  Cc  the  sound  speed)  and  the  outer  radius  a 
is  governed  by 

ma  *  kI2Ct)/a  +  2m  t  Cl  +  T  CT) )  nT, 

with  I  the  current,  Z  the  degree  of  ionization,  n  the  ion  density,  and  T  the 
temperature.  This  is  the  small -skin-depth  limit  of  equation  4  in  Ref.  3, 
i.e.  the  same  as  the  SPLAT  code. 

The  SQUEEZE  code,  originally  written  with  built-in  simple  equations 
for  Aluminum  ionization,  has  been  generalized,  and  takes  over  the  modeling  of 
dynamics  once  the  inner  radius  of  the  annulus  has  shrunk  to  zero.  This  code 
consists  of  solving  simultaneous  ordinary  differential  equations  for  the  radius 
and  temperature  of  a  radiating  cylinder  compressed  by  the  force  of  a  given 
current  I Ct ) . 

Splitting  the  second-order  force  equation  for  r  into  two  first-order 
equations,  we  have: 

r  *  v 

v  »  [A  ZjCT)  T  -  3]  r'1 

T  =  {[ACZ^T)  T  +  BD]  v/r  -  R(r,  T)/r2  +  y/Z^T) 
where  A,  B,  C,  D  and  are  constants  in  the  simplest  model  (the  current  is 
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held  constant  in  that  model),  R(r,  T)  describes  the  temperature-dependence 
and  apacity  dependence  of  the  radiative  cooling.  In  the  optically-thin  limit 
R  is  a  function  of  T  only. 

The  functions  (T )  and  Z^CT)  are  related  to  the  degree  of  ioniza- 

tion  Z(T)  by  Zj(T)  *  +  ZtT)],  describing  the  total  number  density  of 

particles  (ions  and  Z  electrons  per  ion)  with  energy  j  kT,  and  Z^T)  = 

3e./3T  +  j  [1  +  7(T)  +  T3Z/3T],  relating  to  the  partition  of  new  energy  into 
increasing  temperature  and  ionization.  e.(T)  is  the  average  potential  energy 
invested  per  ion  in  stripping  them  to  the  degree  impled  by  electron  temperature 
T. 

The  v  equation  comes  simply  from 

2 

m  r  =  2  it  r  £  (P  -  |^)  with  B  «  ~  (cgs  units) 

N. 

and  with  P  *  (n  +  n.)  T  *  — i-  [1  +  Z(T) ]  T. 
e  1  nr  Z 

The  f  equation  comes  from 


Bt  {l  Ni(1+Z)T  +  Vi}  =  ‘  Ht  {Jm2)  +  (P'b2/8u)  It  {lTr2l)  '  Prad 


+  I-E 


where  the  mean  kinetic  energy 


7^0  1/2  (r' )  2  tt  r'dr'  =  j  m  vZ(r) 

thermalizes  at  a  rate  r  r  with  r  given  from  the  equation  of  motion,  which 
makes  this  term  just  half  the  (P  -  Bj/SiOdV  work.  The  coefficient  AC  thus 


_ j 


24 


relates  to  work  done  by  thermal  pressure  and  BD  relates  to  work  done  by 
magnetic  pressure.  The  Ohmic  heating  term  Is  proportional  to  the  current 
(excluding  any  carried  by  tenuous  outer  plasma)  and  voltage  (.including  inductive 
correction).  The  radiative  cooling  function  R(.r,  t)  is  related  to  the  total 
radiated  power  Prad  by 


R(r.  T)  ( r.  T). 

mr  x  N.  rad 


In  the  optically  thin  limit  this  is  a  function  of  T  only;  in  the  black-body 

3  4 

limit  it  is  proportional  to  r  T  .  The  more  general  form  of  R  in  the  code 
is  to  be  got  by  fits  to  the  MULTICELL  computational  results.  In  spite  of  the 
assumption  of  uniform  density  across  the  cylinder  for  dynamical  purposes,  it 
is  probably  a  better  approximation  in  the  energetics  to  use  radiative  cooling 
results  for  Gaussian-density  cylinders.  Likewise  the  functions  Z  and  3e../3T 
are  calculated  from  smooth-interpolator  fits  to  computational  results,  and  are 
the  same  functions  as  used  in  the  more  sophisticated  codes. ^ 

The  set  of  equations  is  easily  generalized  by  adding  a  circuit 
equation  to  describe  the  current  in  terms  of  the  voltage  and  the  time-varying 
inductance  and  impedance  of  the  plasma,  but  the  resistive  impedance  depends 
on  details  of  the  field  penetration  which  are  not  known  a  priori.  Since 
during  the  short  assembled  "squeeze"  phase  the  current  does  not  change  by  a 
large  percentage,  we  have  chosen  to  take  it  as  an  empirical  constant. 


3.  Variational  Approach  to  Optimizing  Yield 

The  problem  of  optimizing  the  time-integrated  radiation  above  some 
critical  photon  energy  Ec  is,  like  other  optimization  problems,  basically  a 
calcul us-of-variations  problem.  In  the  simple  model  of  a  uniform  or  a 
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Gauss lan-proflle  cylinder  radiator,  there  are,  for  any  given  material,  four 
"coordinates"  satisfying  coupled  first-order  ordinary  differential  equations: 

y^(t)  •  r,  the  cylinder  radius 
yg(t)  »  v,  the  surface  velocity 
y-j(t )  *  T,  the  central  temperature 
y4(t)  *  I,  the  current 

In  the  energy-balance  equation  for  T  (.the  plasma  is  assumed  dense  enough  that 
T  *  T.)  the  radiated  power,  P>  and  P<,  above  and  below  photon  energy  E  , 

Cl  v 

appears  as  a  pair  of  functions  of  n,  r,  and  T.  Conservation  of  mass  allows 
us  to  write  n(r)  trivially  from  initial  conditions,  and  thus  the  power  radiated 
can  be  expressed  as  the  known  functional  forms  P<  (r,  T)  and  P>  (r,  T). 

The  object  of  the  variational  calculation  is  to  maximize 

Y>  =  And  dt  p>  ^ 
zl 

subject  to  the  nonholonomic  constraints  on  r,  T  imposed  by  the  "constraint 
equations" 

r  *  fj(v)  *  v 
v  3  f2(r,  v,  T,  I) 

T  3  f3Cr,  v,  T,  I) 

I  3  f4(r,  v,  T,  I,  V) 

i.e.  yn  3  ^n(yr--y4.  y5).  for  n  *  1...4.  Here  the  independent  "coordinate" 
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y$  *  V(t)  and  the  Initial  conditions  are  varied  In  order  to  achieve  the 
maximum  Y.  V(t)  represents  the  applied  voltage  history. 

The  functions  P<  (r,  T)  and  P>  (r,  T)  correspond  to  functions 
suggested  for  computational  evaluation  In  Section  1.  It  is  expected  that  they 
can  be  fit  by  relatively  simple  algebraic  forms.  The  Integral  Y*  Is  maxi¬ 
mized  by  Introducing  4  Lagrange  multiplier  functions  Xn(t)  and  solving  the 
eight  differential  equations 


X 


n 


4 

Z 

m*l 


3f 
V  — L 
m  3y, 


m 


The  control  function  yg(t)  is  then  determined  by  the  algebraic 


relation 
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X 


+ 


0. 


There  is  the  complexity  of  a  change  of  constraints  at  the  time  tm 
when  the  wires  merge  into  an  annulus,  and  at  the  time  t^  when  the  annulus 
assembles  into  a  cylinder.  One  may  specialize  this  variational  technique  to 
the  assembled  cylinder,  where  the  trajectory  equations  are  those  of  the 
generalized  SQUEEZE  code,  with  the  inclusion,  hopefully,  of  realistic  functions 
P*  (r,  t)  and  P<(r,  T)  fit  from  the  proposed  computer  evaluation  for  Gaussian 
cylinders.  Here,  to  keep  the  system  of  equations  simple,  we  may  to  a  first 
approximation  eliminate  the  circuit  equation  y4  =  f4  (y....yg)  after  the  time 
t, ,  observing  that  the  current  does  not  change  by  a  large  percentage  during  the 
short-lived  assembled  phase  of  the  discharge,  where  the  radiation  is  important. 
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But  one  must  then  couple  this  t  >  t^  dynamics  to  the  t  <  tj  dynamics. 

The  variational  problem,  with  the  lower  Integration  limit  on 
/  P*dt  set  equal  to  tj,  then  formally  provides  an  optimization  of  this 
quantity  by  varying  the  current  and  the  "Initial"  conditions  r^,  v^,  at 
the  beginning  of  the  assembly  (tj).  Since  these  ordinary  differential  equations 
can  be  solved  computationally  much  more  easily  than  the  MHO  equations,  and 
since  they  are  designed  to  lead  to  optimization  equations  rather  than  dynamic 
information,  it  is  hoped  that  they  can  play  a  valuable  role  In  system  design 
at  low  cost,  in  spite  of  their  very  approximate  nature. 


I 
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4.  Caveats  to  These  Methods 


The  radiative  output  of  a  Gausslan-denslty  cylinder  Is  partially 
controlled  by  the  temperature  of  the  outer,  cooler  blanket.  The  radial 
energy-flow  equation  must  contain  a  radl ally-dependent  heat  source  to  give  a 
quasi-steady  radiation  rate  and  Tg  profile.  To  the  extent  self-similar 
compresslonal  heating  provides  the  heat  source,  one  can  take  the  source  to  be 
uniform, l.e.  posltlon-lndpendent.  But  Ohmic  heating,  which  Is  expected  to 
supply  a  smaller  fraction  of  the  energy,  occurs  primarily  near  the  edges,  and 
to  the  extent  that  It  Is  important, It  heats  predominantly  the  blanket  region, 
reducing  its  opacity.  The  radiation  outputs  are  thus  shifted  toward  those  for 
the  isothermal  profile  at  the  same  central  temperature.  This  puts  the  actual 
radiation  picture  somewhere  between  the  isothermal  and  the  blanketed  results 
calculated.  The  sensitivity  to  such  edge  heating  can  be  easily  evaluated. 

Also,  the  heating  and  radiation  are  not  in  balance  for  a  real  time- 
dependent  system.  The  assumption  of  approximate  balance  to  give  a  radiation 
output  above  and  below  photon  energy  e  Is  justified  only  insofar  as  the 

V 

temperature  profiles  are  similar  in  the  dynamic  and  the  steady  problems. 

The  caveats  of  the  combined  WIRES/SQUEEZE  calculation  are  primarily 
in  the  model  equations  and  in  the  zero-dimensional  nature  (i.e.  radial 
structure,  including  time-varying  blanketing  effect, is  not  resolved).  One 
can  see  the  latter  caveat  as  follows:  imagine  two  different  systems  of 
collapsing  annular  plasmas,  one  with  a  sharp  leading  edge  (i.e.  steep  density 
gradient  on  the  inner  edge)  and  the  other  with  a  soft  leading  edge  and  a  steep 
trailing  edge.  When  the  first  plasma  assembles  its  kinetic  energy  is 
'thermalized'  more  rapidly,  pushing  the  radiation  tradeoff  toward  a  shorter 
pulse  at  higher  temperature,  while  the  slower  implosion  may  radiate  longer, 
remaining  cooler.  A  spatially  resolved  MHO  code  can  discriminate  between 
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these  implosions,  while  the  more  approximate  treatment  has  no  consistent  or 
reliable  way  to  do  so.  Nonetheless,  Its  low  cost  makes  It  potentially  quite 
useful  as  a  device  for  obtaining  scaling  with  those  variables  which  It  can 
quantify. 

It  is  simple  to  replace  the  WIRES  run-in  subroutine  with  one 
describing  an  annular  gas  puff.  It  Is  also  simple  to  replace  the  run-in 
radiation/pressure-balance  heating  module  used  in  WIRES  with  an  alternative 
description  of  the  thickness /temperature  relation  for  the  annular  plasma. 

The  obvious  question  is  the  sensitivty  of  results  to  these  models,  which  is 
unknown  as  yet  but  can  be  evaluated  from  runs  with  the  alternate  model  equations. 

Finally,  the  caveats  of  the  calculus-of-variations  approach  to  voltage 
pulse  shaping  and  choice  of  initial  conditions  are  that  such  sets  of  equations 
often  do  not  have  unique  solutions,  or  have  no  reasonable  solutions  without 
further  constraint  on  the  properties  of  the  voltage  pulse.  The  absence  of 
uniqueness,  if  it  occurs,  is  not  serious.  If  there  are  several  equally  good 
solutions  one  simply  picks  the  one  that  is  technologically  easiest.  The  absence 
of  reasonable  solutions,  if  this  occurs,  can  force  one  to  complicate  the 
calculation  by  add i tonal  constraints. 

The  variational  approach,  because  it  uses  the  models  of  calculations 
(1)  and  (2),  also  suffers  from  their  potential  drawbacks;  nonetheless,  it  is 
a  significant  step  forward  at  relatively  low  cost,  and  we  do  not  expect  these 
drawbacks  to  be  determining  factors  in  system  design. 
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C.  Some  Recommendations  for  Experiments 

During  1981  JAYCOR  also  suggested  experiments,  at  the  request  of 
NRl,  which  could  be  done  on  present  devices  and  which  might  explore  techniques 
for  producing  energetic  radiation  nonthermally  where  the  plasma  bulk  tempera¬ 
ture  cannot  be  raised  enough  to  emit  such  radiation  thermally. 

Two  experiments  were  proposed  on  diode-imploded  hlgh-z  plasma  loads. 
The  first  Involves  Implosion  of  an  array  of  low-mass  silicon  wires  onto  an 
aluminum  wire,  and  Implosion  of  an  outer  array  of  silicon  wires  with  an  inner 
array  of  aluminum.  The  objectives  of  this  experiment  would  be  to  observe 
radiation  or  any  enhancement  by  the  atomic  and  radiation  physics  specific  to 
these  two  materials,^®  to  compare  the  radiation  with  and  without  a  central 
wire, ^  and  to  use  the  observed  implosion  time  and  diode  current  and  voltage 
traces  to  infer  the  degree  of  magnetic  shielding  of  inner  conductors  by  outer 
ones,  on  the  voltage  pulse  timescale. 

The  second  experiment  involves  varying  the  charging  voltage  applied 
to  a  conventional  low-mass  array  of  high-z  wires  to  observe  and  study  the 
threshold-like  onset  of  hard  x-ray  (10-100  keV)  radiation  from  small  dense 
regions  of  the  final  pinched  discharge,  as  observed  by  x-ray  pinhole  photo¬ 
graphy.^  A  first  objective  here  would  be  to  correlate  the  voltage  at  which 
this  threshold  occurs  with  the  theory  of  modified-two-stream  unstable  regions 

IQ 

of  low  density  in  series  with  higher  density  beads.  The  theory  predicts 
that  if  the  sausage  instability  reduces  the  density  in  the  necked-off  regions 
sufficiently,  enhanced  hard  radiation  will  be  produced  as  electrons  heated  in 
the  low-density  "beads".  A  second  objective  of  this  experiment  is  to  determine 
whether  this  phenomenon  can  be  used,  as  predicted,  to  tailor  the  "hard" 
component  of  the  spectrum  and  make  significant  radiation  production  in  the 
10-100  keV  range  from  discharges  too  massive  or  too  high  in  Z  to  reach  the 

temperatures  needed  for  significant  thermal  K-line  radiation. 
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1. 


ualltative  Physics  of  Double  Arrays 

Implosions  of  an  annular  array  of  silicon  fibers  onto  an  aluminum 
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center  wire,  or  of  a  large-radius  silicon  array  and  a  smaller  aluminum  array 
on  the  axis,  appear  Interesting.  The  atomic  and  radiative  Interactions  of 
the  silicon  and  aluminum  at  high  temperatures  have  been  treated  In  a  1981 
report.^  Considerations  here  are  that  the  temperatures  of  the  SI  and  A1 
components  must  become  more  or  less  comparable  at  the  same  time.  Although 
physical  proximity  of  the  two  materials  Is  not  explicitly  required,  the 
requisite  temperatures  (MOO  eV)  are  created  only  during  the  peak  compression 
on  axis,  so  that  for  practical  purposes  the  more  interesting  and  efficient 
radiation  process  can  come  about  only  when  the  load  collapses  onto  the  axial 
wire,  or  when  the  two  loads  coincide  on  axis.  Separate  arrays  of  different 
mass  will  tend  to  implode  at  different  rates,  corrected  somewhat  by  magnetic 
shielding  of  the  inner  load  by  the  outer  one,  so  that  the  implosion  timescale 
will  tend  to  be  that  of  the  slower  load. 

While  an  annular  outer  conductor  shields  the  applied  electromagnetic 
field  from  an  array  of  conductors  inside  it,  this  is  not  necessarily  true  if 
the  outer  conductor  is  a  small  number  of  wires  distributed  in  annular  fashion. 
To  the  extent  that  the  elctromagnetic  pulse  is  not  transmitted  to  the  inner 
array,  the  time-history  of  such  double  array  implosion  with  staggered  wires 
would  be  expected  to  look  like  Fig.  1,  in  which  the  outer  array  of  implodes 
past  the  inner  one  but  loses  its  J  x  8  driving  force  in  so  doing,  while  the 
second  array  accelerates  from  rest.  By  suitable  choice  of  the  array  masses 
and  radii,  the  two  implosions  can  be  made  to  coincide  on  reaching  the  axis. 

When  the  magnetic  shielding  is  poor,  the  inner  array  is  subject  to 
some  J  x  B  force  before  being  overtaken,  and  the  timing  of  the  two  implosions 
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Is  altered.  This  provides  a  potentially  sensitive  diagnostic  for  the  distribu¬ 
tion  of  J  x  B  forces  In  such  arrays,  and  should  shed  light  on  the  degree  of 
electromagnetic  shielding  (and  hence  the  coupling  efficiency)  in  the  conventional 
annular  wire  arrays. 

2.  Observations  on  Beaded  Discharges 

There  is  a  reasonable  amount  of  evidence  showing  that  at  high 
voltages,  exploding  wires  of  high  atomic  number  exhibit  beading,  including 
millimeter-size  blobs  of  dense  plasma  as  seen  in  the  thermal  radiation  spectrum, 
and  smaller  hot  spots  of  energetic  x-ray  emission.  These  discharges  are  un¬ 
stable  to  sausage  instability,  and  it  is  reasonable  to  expect  the  observed 
behavior.  Imploded  annular  arrays  of  wires  exhibit  a  more  complex  variety 
of  inhomogeneities,  including  beading  and  helical  distortions,  reflecting 

the  fact  that  they  are  unstable  to  corkscrew-like  modes  in  addition  to  the 
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m=0  sausage.  It  has  been  proposed  that  these  inhomogeneities  can  be  useful 
in  allowing  intense  radiation  by  heating  electrons  strongly  in  the  low-density 
regions  and  bringing  about  the  energy  deposition  at  the  surface,  or  in  the 
volume,  of  the  denser  beads. 

While  the  beading  is  a  quite  chaotic  process,  certain  conditions 

on  voltage  and  minimum  density  on  axis  are  predicted  to  give  rise  to  such 

strong  electron  heating.  This  phenomenon  has  not  yet  been  explored  experimen- 

18 

tally  except  in  single  wire  loads  and  then  no  attempt  was  made  to  study  the 
scaling  of  the  hard  radiation  with  voltage  or  to  estimate  densities  in  the  less 
dense  regions. 
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3.  Potential  Systems  Impact  of  the  Experiments 

The  simulation  community  Is  presently  not  able  to  build  a  device 
with  the  required  combination  of  total  radiated  power,  spectrum,  and  timescale 
for  DNA  objectives.  The  importance  of  making  hlgh-z  diode-imploded  plasmas 
radiate  nonthermally  in  theK-line  and  Bremsstrahlung  continuum  is  that 
insufficient  energy  Is  available  with  present-generation  devices  to  make  them 
radiate  in  this  spectral  range  thermally,  unless  unreasonably  low-mass  arrays 
are  envisioned.  Two  nonthermal  processes,  discussed  in  reports  but  not  yet 
explored  experimentally,  hold  some  hope  for  enhancing  the  radiation  from 
plasma  sources  to  suitable  levels.  The  experiments  proposed  were  designed  as 
preliminary  looks  at  those  two  processes,  and  in  passing  provide  important 
information  about  the  electrodynamics  and  atomic  physics  of  ongoing  experiments. 
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